Unwinding dynamics of double-stranded polymers 
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" Abstract 

We consider the unwinding of two lattice polymer strands of length N which are initially wound 
around each other in a double helical conformation and evolve through Rouse dynamics. The 
problem relates to quickly bringing a double-stranded polymer well above its melting temperature, 

-a ■ 

i.e., binding interactions between the strands are neglected, and the strands separate from each 

O : 

other as it is entropically favorable for them to do so. The strands unwind by rotating around each 
7-H ■ other until they separate. We find that the process proceeds from the ends inwards; intermediate 

>: 

conformations can be characterized by a tightly wound inner part, from which loose strands are 

m . 

sticking out, with length I ~ i - 39 . The total time needed for the two strands to unwind scales as 

^ ; 

a power of N as t u ~ jy2.57±o.03^ We present a theoretical argument which suggests that during 

! this unwinding process, these loose strands are far out of equilibrium. 

t— I . 

> '■ 

.!-( , PACS numbers: 36.20.Ey, 87.15.H-, 87.14.gk 
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I. INTRODUCTION 



There are several known examples of polymers in nature that are composed by two or more 
strands arranged in a helical conformation. One is the double helical structure of DNA, in 
which two complementary sequences are held together by hydrogen bonding between A/T 
or C/G nucleotides. Another example is the triple helical structure of collagen. Under 
the appropriate thermodynamic conditions the multi-stranded structure becomes unstable 
and the strands dissociate from each other. For a DNA molecule, the dissociation of the 
two strands is usually referred to as the melting transition or denaturation; this happens 
when the temperature is increased typically above 80° C in standard conditions pj. In 
view of its importance in many biological and biotechnological processes, various models of 
DNA melting have been developed to study the thermodynamical behavior as function of 
the sequence .ength and composition. Mode.s originaHy proposed by Po,and and Scheraga 
2] (PS) and by Peyrard and Bishop [3] (PB) take into account the different unbinding 
energies required to dissociate complementary G/C or A/T pairs using various types of 
approximations, yet they are simple enough to deal with very long DNA sequences. Using 
an appropriate choice of parameters they both fit reasonably well experimental data for the 
melting temperature of DNA sequences {4]. 



The melting dynamics has also been studied. In 1986, Baumgartner and Muthukumar 
performed Monte Carlo simulations of the disentangling of two initially intertwined chains. 
After an initial softening of the original double helix, they reported that the time required 
for the actual unraveling of the chains scales with polymer length as a power law with 
exponent 3.3 ± 0.2. Computer technology at that time allowed for simulations of 2, 4 and 8 
turns only, and their fitted exponent might very well change with increasing chain length. 
The dynamics of the opening of bubbles in an entangled chain has also been studied (]- 
[^l. In addition different extensions of the PS and PB models were considered to include 
helical degrees of freedom of DNA molecules 10Ml3[] . These extensions allow for rapid 
computations of equilibrium and dynamical properties of the melting transition, at the cost 
of introducing approximations which are poorly controlled and sometimes even questionable. 
The moves introduced to update PS models usually neglect the helicity of dsDNA, which is 
locally conserved due to impenetrability of the two strands (the unbinding of the two strands 
forming the double helix requires a release of the twist through a rotation of these strands 
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with respect to each other). With the denaturation times characterized by the scaling law 
Td ~ N 13 for the DNA strand length N (in base pairs), the resulting values of /3 thus range 
from /3 = js| to ft = 4/3 Q. If local moves preserving the DNA helicity are introduced, a 
very slow melting is instead observed, with /5 ~ 3 [9| . A second simplification intrinsic in PS 
models is that helical fragments and loops are described by equilibrium partition functions. 
This description can be too simplified in systems where the dynamics of the unwinding 
process is too fast to allow for the full structural relaxation within these loops. PB models 
have as a major simplification the possibility for complementary bases to orbit around a 
virtual central axis of the macromolecule. This yields denaturation times scaling linearly 
with the chain length [3]. Given such spread of results, and uncontrolled simplifications in 
these models, we must conclude that we are still missing the understanding of how DNA 
length affects the time for its strands to disentangle. 

In order to gain more direct insight into the dynamics of the melting process of 
double-stranded molecules we investigate the unwinding dynamics of double-stranded three- 
dimensional long polymers using Monte Carlo simulations in the absence of hydrodynamic 
interactions. No binding energies between the two strands are taken into account during the 
unwinding process, corresponding to the case of a double helix brought rapidly to a temper- 
ature well above its melting point. For such a setup, the strands unwind from each other as 
it is entropically favorable for them to do so. We follow a procedure very similar to that of 
Ref. 5|, except that our chains are much longer, and more tightly wound. Figure [T] shows 
three configurations: (a) at the early stages of unwinding, (b) during the unwinding process, 
and (c) at the end of the unwinding, when the two stands are separating from each other. 
The main scope of this manuscript is two-fold: in section [III we study the unwinding time 
t u as a function of iV — high precision simulations for polymers of length up to N = 1000 
show that t u scales as a power-law t u ~ N 13 with j3 = 2.57 ± 0.03. In section [TTT1 we show 
that the intermediate conformations can be characterized as a tightly wound inner part, to 
which unwound single strands are connected. The length of the unwound strands increases 
with time in a powerlaw fashion as l(t) ~ t 39 . With a theoretical argument, we find an 
upper bound for the radius of gyration of the loose ends, which excludes the equilibrium 
value. This demonstrates that the unwinding is a far from equilibrium process. 
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II. UNWINDING TIME 



In the simulation the polymers reside on a face-centered-cubic lattice with a lattice spacing 
of \/2, and are initialized in a double-helical state. The polymers evolve in time through a 
long sequence of single-monomer moves, under the restriction that at all times, the polymer 
backbones are self- and mutually avoiding. Each allowed move occurs with a statistical 
rate of unity. To give the polymers some elasticity, the self-avoidance condition is lifted for 
monomers which are direct neighbors along the same chain. A detailed description of this 
lattice polymer model, its computationally efficient implementation, and a study of some of 
its properties and applications can be found in [l4]. This model reproduces known features 
of the Rouse dynamics 15[ and of the equilibrium properties [16] of single self-avoiding 
polymers. As the moves respect the no-crossing condition between strands, we expect that 
the long time behavior of unwinding discussed in this study is of universal nature and is not 
affected by microscopic details and lattice effects. 

Let rf\t) and r^ 2 \t) be the lattice positions of the z-th monomers on the two strands 
at time t (0 < i < N). We consider the minimal distance between two strands defined as 
d m i n (t) = min r^\t) — rf\t) . The inset in Figure [2] shows a plot of d^ ain (t) as a function of 
time for a run. The choice of an initial double-helical conformation implies that d m i n = y2 
at t — 0. Note that d m i n (t) remains constant up to a time t ~ 1.9 ■ 10 6 in the inset of Fig. [2] 
and then starts fluctuating and increasing in time. We define the unwinding time t u as the 
time at which d m i Q (t) exceeds some threshold value for the first time. For the threshold value 
do we took d$ = 10 and d^ = 20. The higher threshold value gives a slightly higher estimate 
of the unwinding time [r^ 10 '' and r^ 20 ^ in Fig. [2] . However, as the polymer length increases 
the ratio of r^ 10 ^ /r« 20 ^ converges to 1 as shown in Table [H hence the two quantities have the 
same scaling behavior in N. 

Figure [2] and Table [T] show the behavior of ri 10 '' and function of polymer length. 

We note that the scaling of unwinding times is a power of the strand length: From a 
linear regression of the data for iV > 30 we find the values /3 = 2.58 ± 0.03 (d = 10) and 
(3 = 2.56±0.03 (d = 20), from which we obtain the result anticipated above j3 = 2.57±0.03. 
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III. CHARACTERIZATION OF INTERMEDIATE CONFORMATIONS 



Given the topological constraint each strand faces while unwinding, we expect the unwind- 
ing dynamics to unroll from the two ends of the initial double-stranded complex, progressing 
inwards with increasing time (Fig. [3]). Note that because of the elasticity of the model used 
in the simulations, a partial opening up of the inner wound part is not ruled out by the 
model: it is the physics of the problem that seems to suppress this. In order to connect 
this physical picture with the observed scaling t u ~ ]y2.57±o.o3 we considered the quantity 
dm_in{ht) — m hi ^1 if) ~ "if if) which is the minimal distance from the z-th monomer of 
the first strand to any other monomer of the second strand. Its average square (d^ in (i, £)) is 



plotted in Fig. H] for a strand length equal to N = 500. The different data are for increasing 
time step snapshots taken at time intervals equal to t/10 5 = 1, 2, . . . , 10 (from bottom to 
top). The quantity (d^^ji, t)) is minimal in the middle, while it increases in time from the 
two edges, in agreement with the physical picture proposed in Fig. [31 We consider now the 
normalized profile (^min(M))/ (^min(O)O)- For this quantity we expect the following scaling 
behavior, as a function of the distance from the end monomer (i = 0) 

with /() a scaling function and l(t) a characteristic length depending on time. The inset of 
Fig. H]show that the normalized profiles at different times collapse when a rescaling i/t 0,39 
is used, which implies I ~ t ' 39 . This is consistent with the exponent determined from the 
scaling of unwinding time as I ~ t 1 ^ = t a39 . 

In order to gain insight into the unwinding we set up a simple analytical model of the 
process, assuming that the unwinding is sufficiently slow so that the conformation of the 
loose strands can be approximated by equilibrium ones at all times. We will show that 
this approach predicts an unwinding dynamics which is slower than what is observed in 
simulations. We therefore conclude that the unwinding we observe is a far from equilibrium 
process. 

Consider an intermediate conformation which consists of two single-stranded coils of 
N — s(t) attached to a double stranded helical part of length s(t). We estimate the free 
energy F as a function of s(t) from the partition function for a double helical segment of 
length s is Zh e ii x ~ /J>% and that of a single stranded coil Z co n ~ ^~ s , where \i are the 
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connectivity constants. In other words 

0F = -s]ogQM h )+2(N-s)]og(fM c ). (2) 

As an infinitesimal portion ds of the double helix unwinds, the change in free energy is thus 
given by 

dF = —f3^ 1 \og(Hc / Hhjds = —Kids. (3) 

During this infinitesimal unwinding, the single-stranded coils are displaced over a distance 
dr due to the rotational motion around the axis of the helix: the coils describe a fraction of 
a circle of radius R v perpendicular to the axis of the helix, where _R„ is the distance from 
the helix axis to the coil's center of mass. In equilibrium we expect R v ~ (N — s) u , with v 
the Flory exponent. This implies that 

dr oc R v ds oc (N — s) u ds. (4) 

During this process the work done against friction equals 

dW = jrdr = ^R 2 v s ds, (5) 

where for Rouse dynamics the friction on the single-stranded coils is proportional to their 
lengths: 7 oc iV — s. Since the work done against friction cannot exceed the available free 
energy, we obtain the inequality 

K 2 (N - s) 1+2u s ds < Kxda. (6) 

If this inequality were saturated, the unwinding process would be described by 

K./K^-lN-sit)}^^. (7) 

The integration of Eq. (j7J) yields a scaling of the unwinding time as r u ~ j\r 2l/+2 = N 3AS 
(since v = 0.59), which is obtained from Eq. ([7]) by setting s(r u ) = N. Indeed, a visual 
investigation of the unwound parts during simulations reveals that these show a tendency 
to be spiral-like, and contain a significant amount of coiling. A precise quantification of the 
amount of spiraling is very difficult due to the lack of a clean definition of the central axis of 
the spiral. This argument predicts a very slow unwinding compared to that which is actually 
observed in simulations, suggesting that unwinding proceeds through non-equilibrium states. 
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In an early study of unwinding [5] the scaling r ~ jV 3 - 3 ^ 2 ) was computed for shorter polymers 
(up to iV = 65; with 2, 4 or 8 turns in the double helix) and for helices less tightly bound than 
considered here. The exponent of [5| most likely describes a pre-asymptotic scaling regime. 
The exponent reported in {$]] is however compatible with that predicted from Eq. (J7|). This 
suggests that the early stages of unwinding are probably well captured by the mechanisms 
underlying the derivation of Eq. (jZ]), but that the asymptotic scaling regime is dominated 
by a mechanism which is faster for longer chains. 



IV. SUMMARY AND CONCLUSIONS 



Summarizing, we introduced a lattice model for studying the unwinding dynamics of a 
long three-dimensional double stranded polymer, with excluded volume effects taken into 
account. The lattice nature of the model, combined with an efficient encoding of the dynam- 
ics, allows one to simulate long polymers (up to N = 10 3 ) for very long time (t — 10 7 Monte 
Carlo steps). Our numerical results show that the unwinding time scales with the polymer 
length as a power-law with exponent /3 = 2.57(3) with Rouse dynamics. An analysis of a 
simple analytical model of the process suggests that the unwinding we observe is a far from 
equilibrium process, therefore it cannot be understood in terms of a slow dynamics evolving 
through quasi-equilibrium states, e.g. using free energy arguments. 
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N 


T (io) 


T (20) 


r (20) / (10) 


60 


1.23 


10 4 


1.12 


10 4 


1.093 


80 


1.80 


10 4 


1.97 


10 4 


1.096 


100 


3.04 


10 4 


3.29 


10 4 


1.082 


120 


4.73 


10 4 


5.06 


10 4 


1.070 


150 


8.27 


10 4 


8.75 


10 4 


1.055 


200 


1.81 


10 5 


1.90 


10 5 


1.048 


300 


4.97 


10 5 


5.16 


10 5 


1.038 


400 


1.03 


10 6 


1.06 


10 6 


1.032 


500 


1.82 


10 6 


1.88 


10 6 


1.036 


600 


3.00 


10 6 


3.10 


10 6 


1.032 


800 


6.22 


10 6 


6.35 


10 6 


1.021 


1000 


1.11 


10 7 


1.13 


10 7 


1.017 



TABLE I: Average unwinding times r« and function of polymer length N for threshold 

values cIq = 10 and cIq = 20, and the ratio of these two times. Times are obtained by averaging 
over 120 simulations. 
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FIG. 1: (Color online) Snapshots of the polymer configurations during unwinding for two strands 
of length N = 100 each. The initial conformation is fully double helical all along its length, (a) 
Snapshot after short time from the beginning of the simulation; opening begins mainly from the 
two ends, although small bubbles within the chain are also visible, (b) Snapshot at later times, 
(c) Separation. 
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2 3 

10 10 

N 

FIG. 2: Double-logarithmic plot of the average unwinding time as a function of strand length. 
The circles are obtained with an unwinding threshold of c?q = 10, while the crosses correspond to 
a threshold of = 20. The straight dashed line is a fit to the data corresponding to an unwinding 
exponent of /3 = 2.57. Inset: Plot of (^^(t) vs. time for a run with strands of length N = 500. 
The arrows indicate the first time that this distance reaches its threshold value d% = 10 or 20. 




FIG. 3: Sketch of a double-stranded polymer during the unwinding dynamics. At time t we expect 
to find a double stranded region of curvilinear length s(t) terminates with two single strands of 
lengths (N — s{t))/2 at both edges, s(t) being a decreasing function of t. 
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i 

FIG. 4: Average squared minimal distance {df nin (i,t)) as a function of monomer number i, for 
polymers with length N = 500. From bottom to top, the curves are obtained at times i/10 5 = 
1,2... 10. Inset: Collapse of (t4in(M))/(<4i n ( >*)) for polymers of length N = 1000 again at 
times t/W 5 = 1, 2 ... 10. The vertical scale is normalized so that all curves start from a common 
value at i = 0. The horizontal scale is divided by t a39 . This exponent is consistent with that 
obtained from the scaling of the unwinding time. 
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